#Datos sobre pobreza en los municipios de Mexico: Indice de rezago educativo, pobreza, indice de carencia en acceso a servios de salud
#y indice de carencia en alimentos
# data table
import pandas as pd
import os

datamuni=pd.read_excel(os.path.join('data','indicadores de pobreza municipal, 2015 (3).xlsx'),
                     dtype={'Ubigeo': object})
datamuni.head()

datamuni.info()

#Leemos el mapa que contiene los municipios de Mexico
import os, geopandas as gpd

MunicipiosMex=gpd.read_file(os.path.join("maps","mapa_mexico","Division_Municipal_Mexico_2010.shp"))
MunicipiosMex
                            

datamuni

MunicipiosMex

#Calcular la cantidad de filas y de columnas 
len(datamuni.Municipio),len(MunicipiosMex.NOM_MUN)

#Verificamos si hay repetidos
datamuni.Municipio.duplicated().sum(),MunicipiosMex.NOM_MUN.duplicated().sum()

#La nueva columna se ve así:
datamuni['estaMuni'].head()

#Esto es lo que no se pudo igualar:
len(nomatch_df), len(nomatch_gdf)

# elija la coincidencia más cercana de nomatch_gdf para un valor en nomatch_df 
from thefuzz import process
[(dis,process.extractOne(dis,nomatch_gdf)) for dis in sorted(nomatch_df)]

#Como calculo el diccionario 
{dis:process.extractOne(dis,nomatch_gdf)[0] for dis in sorted(nomatch_df)}

#Verificamos
MunicipiosMex.info()

#Eliminamos los que no vamos a utilizar
bye=['Municipio', 'Entidad_Federativa','CVE_EDO', 'CVE_MUNI']
MunicipiosMex.drop(columns=bye,inplace=True)

# keeping
MunicipiosMex.head()

#Explorando una variable
#Esta vez, exploramos estadísticamente una variable en el mapa:
#Estadistica: Pobreza
MunicipiosMex.IC_Reza_Educ.describe()

#Graficamos los datos estadistucos, histogramas
import seaborn as sea

sea.displot(MunicipiosMex.IC_Reza_Educ, color='yellow')

#This plot tells you the distribution of the values, but not the presence of outliers, which you are revealed in a boxplot:
#Diagrama de cajas
sea.boxplot(MunicipiosMex.IC_Reza_Educ, color='yellow',orient='h')

# a frequency table by default/ es un algoritmo que tu no le dices cuantos intervalos quieres, y evalua la similutd 
ht

# group label
ht.yb

# labels and counts
np.unique(ht.yb,return_counts=True)

# ht.yb into a pandas Series convertilo a pandas 

pd.Series(ht.yb).value_counts()

# these are the cuts, but it is not including the min value/ El valor donde hago el corte de los valores. 
#No dice el minimo valor, pero si el maximo
ht.bins

#Graficamos
sea.displot(MunicipiosMex.IC_Reza_Educ, bins=HT_bins,color='yellow')

#Todo lo que hemos calculado, el adcm de todas las tecnicas que hemos calculado y solo nos quedamos con las 3 que tengan menor error
class5 = q5, ei5,msd, ht, mb5, fj5, jc5, mp5
# Collect ADCM for each classifier
fits = np.array([ c.adcm for c in class5])
# Convert ADCM scores to a DataFrame
adcms = pd.DataFrame(fits)
# Add classifier names
adcms['classifier'] = [c.name for c in class5]
# Add column names to the ADCM
adcms.columns = ['ADCM', 'Classifier']
ax = sea.barplot(
    y='Classifier', x='ADCM', data=adcms, palette='Pastel1'
)

#Cuantos distritos hay en cada caso/ haow many cases per scheme:

class5new = q5,fj5, jc5
pd.DataFrame(
    {c.name: c.counts for c in class5new},
    index=['Class-{}'.format(i) for i in range(5)]
)

#Voy a graficar una de ellas.
import matplotlib.pyplot as plt

f, ax = plt.subplots(1, figsize=(9, 9))
MunicipiosMex.plot(column='QU', 
        cmap='viridis', 
        categorical=True,
        edgecolor='white', 
        linewidth=0., 
        alpha=0.75, 
        legend=True,
        legend_kwds=dict(loc=2),
        ax=ax
       )

ax.set_axis_off()

#Voy a graficar la segunda de ellas.
import matplotlib.pyplot as plt

f, ax = plt.subplots(1, figsize=(9, 9))
MunicipiosMex.plot(column='FJ', 
        cmap='plasma', 
        categorical=True,
        edgecolor='white', 
        linewidth=0., 
        alpha=0.75, 
        legend=True,
        legend_kwds=dict(loc=2),
        ax=ax
       )

ax.set_axis_off()

#Voy a graficar una de ellas.
import matplotlib.pyplot as plt

f, ax = plt.subplots(1, figsize=(9, 9))
MunicipiosMex.plot(column='JC', 
        cmap='magma', 
        categorical=True,
        edgecolor='white', 
        linewidth=0., 
        alpha=0.75, 
        legend=True,
        legend_kwds=dict(loc=2),
        ax=ax
       )

ax.set_axis_off()

#insatisfechas
selected_variables = ['IC_Reza_Educ',
                     'IC_Care_en_Acce_Servi_Salud', 
                     'IC_Care_en_Ali']
MunicipiosMex[selected_variables].describe()

sea.boxplot(MunicipiosMex[selected_variables])

#
MunicipiosMex[selected_variables].hist();

#Mirada bivariada. Verificar que las tres variables tengan la misma correlacion o todas posi o todas nega. Todas presentan correlacion positiva
sea.pairplot(
    MunicipiosMex[selected_variables], kind="reg", diag_kind="kde"
)

pd.DataFrame(rs_result).boxplot(vert=False)

pd.DataFrame(pt_result).boxplot(vert=False)

#Elegimos el que tenga menos cantidad de datos atipicos, este será el elegido 
pd.DataFrame(qt_result).boxplot(vert=False)

f, axs = plt.subplots( ncols=3, nrows=2,figsize=(10, 10))
axs = axs.flatten()
# Start a loop over all the variables of interest
for i, col in enumerate(selected_variables_t + selected_variables):
    # select the axis where the map will go
    ax = axs[i]
    # Plot the map
    MunicipiosMex.plot(
        column=col,
        ax=ax,
        scheme="Quantiles",
        linewidth=0,
        cmap="RdPu",
    )
    # Remove axis clutter
    ax.set_axis_off()
    # Set the axis title to the name of variable being plotted
    ax.set_title(col)
# Display the figure
plt.show()

qt_result #Es el dato que voya utilizar con menoserrores 

from scipy.cluster import hierarchy as hc


Z = hc.linkage(qt_result, 'ward')
# calculate full dendrogram
plt.figure(figsize=(25, 10))
plt.title('Hierarchical Clustering Dendrogram')
plt.xlabel('cases')
plt.ylabel('distance')
hc.dendrogram(
    Z,
    leaf_rotation=90.,  # rotates the x axis labels
    leaf_font_size=1,  # font size for the x axis labels
)
plt.show()

MunicipiosMex["hc_3"].value_counts()

#Promedio de grupo por cada variable
MunicipiosMex.groupby("hc_3")[selected_variables].mean()

MunicipiosMex[selected_variables]

#Preparamos una tabla con los datos obtenidos. Transformamos la data en formato long
# Index db on cluster ID
MunicipiosMex_long = MunicipiosMex.set_index("hc_3").copy()
# Keep only variables used for clustering
MunicipiosMex_long = MunicipiosMex_long[selected_variables]
# Stack column names into a column, obtaining
# a "long" version of the dataset
MunicipiosMex_long = MunicipiosMex_long.stack()
# Take indices into proper columns
MunicipiosMex_long = MunicipiosMex_long.reset_index()

# Check out result
MunicipiosMex_long.head()

# Rename column names
MunicipiosMex_long = MunicipiosMex_long.rename(
    columns={'level_1': 'Attribute', 0: 'Values'}
)
MunicipiosMex_long.head()

# Configurar las facetas para poder hacer esta grafica de perfilado
facets = sea.FacetGrid(
    data=MunicipiosMex_long,
    col="hc_3",
    hue="Attribute",
    sharey=False,
    sharex=False,
    aspect=2,
    col_wrap=3,
)
# Construye la trama como una `sns.kdeplot`
facets.map(sea.kdeplot, "Values", fill=True).add_legend()

#Graficamos el mapa 
# Configuramos la figura
f, ax = plt.subplots(1, figsize=(9, 9))
#Trzamos valore unicos
MunicipiosMex.plot(
    column="hc_3", categorical=True, legend=True, linewidth=0, ax=ax
)
# Remove axis
ax.set_axis_off()
# Display the map
plt.show()

# first one. Cuales son sus vecinos 
MunicipiosMex.head(1)

# cantidad de vecinos de ese distrito
w_rook.neighbors[0] 

# see
base=MunicipiosMex[MunicipiosMex.NOMEDO=="BAJA CALIFORNIA"].plot()
MunicipiosMex.iloc[w_rook.neighbors[0] ,].plot(ax=base,facecolor="yellow")
MunicipiosMex.head(1).plot(ax=base,facecolor="red")

w_queen.neighbors[0]

base=MunicipiosMex[MunicipiosMex.NOMEDO=="BAJA CALIFORNIA"].plot()
MunicipiosMex.iloc[w_queen.neighbors[0] ,].plot(ax=base,facecolor="yellow")
MunicipiosMex.head(1).plot(ax=base,facecolor="red")

w_knn.neighbors[0]

base=MunicipiosMex[MunicipiosMex.NOMEDO=="BAJA CALIFORNIA"].plot()
MunicipiosMex.iloc[w_knn.neighbors[0],].plot(ax=base,facecolor="yellow")
MunicipiosMex.head(1).plot(ax=base,facecolor="red")

# all the neighbors by row. Todos los vecinos que tiene cada fila
w_knn.neighbors

# the matrix of neighboorhood: 

pd.DataFrame(*w_knn.full()).astype(int) # 1 means both are neighbors

#cuenta de ceros de cuantos no vecinos hay
w_knn.nonzero

# pct of neighboorhood (density)
w_knn.pct_nonzero

# Hay algun municipio sin vecino?
w_queen.islands

from esda.moran import Moran #sI EL VALOR DE UN DISTRITO ESTA CORRELACIONADO CON SU DISTRITO VECINO

moranICED = Moran(MunicipiosMex['IC_Reza_Educ'], w_knn)
moranICED.I,moranICED.p_sim

# Grafico de la correlcion con el ICED general con el ICED de los vecinos
from splot.esda import moran_scatterplot

fig, ax = moran_scatterplot(moranICED, aspect_equal=True)
ax.set_xlabel('ICED_std')
ax.set_ylabel('SpatialLag_ICED_std')
plt.show()

from esda.moran import Moran #sI EL VALOR DE UN DISTRITO ESTA CORRELACIONADO CON SU DISTRITO VECINO

moranSALUD = Moran(MunicipiosMex['IC_Care_en_Acce_Servi_Salud'], w_knn)
moranSALUD.I,moranSALUD.p_sim

# Grafico de la correlcion con el ICSalud general con el ICSalud de los vecinos
from splot.esda import moran_scatterplot

fig, ax = moran_scatterplot(moranSALUD, aspect_equal=True)
ax.set_xlabel('SALUD_std')
ax.set_ylabel('SpatialLag_SALUD_std')
plt.show()

from esda.moran import Moran #sI EL VALOR DE UN DISTRITO ESTA CORRELACIONADO CON SU DISTRITO VECINO

moranAli = Moran(MunicipiosMex['IC_Care_en_Ali'], w_knn)
moranAli.I,moranAli.p_sim

# Grafico de la correlcion con el ICSalud general con el ICSalud de los vecinos
from splot.esda import moran_scatterplot

fig, ax = moran_scatterplot(moranAli, aspect_equal=True)
ax.set_xlabel('Ali_std')
ax.set_ylabel('SpatialLag_Ali_std')
plt.show()

# El diagrama de dispersión con información local
from esda.moran import Moran_Local

# calculate Moran_Local and plot
lisaICED = Moran_Local(y=MunicipiosMex['IC_Reza_Educ'], w=w_knn,seed=2022)
fig, ax = moran_scatterplot(lisaICED,p=0.05)
ax.set_xlabel('ICED_std')
ax.set_ylabel('SpatialLag_ICED_std')
plt.show()

# el mapa con los puntos y los valores atípicos

from splot.esda import lisa_cluster
f, ax = plt.subplots(1, figsize=(12, 12))
plt.title('Spots and Outliers')
fig = lisa_cluster(lisaICED, 
                   MunicipiosMex,ax=ax,
                   legend_kwds={'loc': 'center left', 
                                'bbox_to_anchor': (0.7, 0.6)})
plt.show()

# significance, la significancia de cada municipio
lisaICED.p_sim

# quadrant: 1 HH,  2 LH,  3 LL,  4 HL
pd.Series(lisaICED.q).value_counts()

#Esta data la metemos al dataFrame
MunicipiosMex['ICED_quadrant']=[l if p <0.05 else 0 for l,p in zip(lisaICED.q,lisaICED.p_sim)  ]
MunicipiosMex['ICED_quadrant'].value_counts()

#Graficamos
labels = [ '0 no_sig', '1 hotSpot', '2 coldOutlier', '3 coldSpot', '4 hotOutlier']

MunicipiosMex['ICED_quadrant_names']=[labels[i] for i in MunicipiosMex['ICED_quadrant']]

MunicipiosMex['ICED_quadrant_names'].value_counts()

from matplotlib import colors
myColMap = colors.ListedColormap([ 'ghostwhite', 'red', 'green', 'black','orange'])



# Set up figure and ax
f, ax = plt.subplots(1, figsize=(12,12))
# Plot unique values choropleth including
# a legend and with no boundary lines

plt.title('Spots and Outliers')

MunicipiosMex.plot(column='ICED_quadrant_names', 
                categorical=True,
                cmap=myColMap,
                linewidth=0.1, 
                edgecolor='white',
                legend=True,
                legend_kwds={'loc': 'center left', 
                             'bbox_to_anchor': (0.7, 0.6)},
                ax=ax)
# Remove axis
ax.set_axis_off()
# Display the map
plt.show()

# The scatterplot with local info
from esda.moran import Moran_Local

# calculate Moran_Local and plot
lisaSALUD = Moran_Local(y=MunicipiosMex['IC_Care_en_Acce_Servi_Salud'], w=w_knn,seed=2022)
fig, ax = moran_scatterplot(lisaSALUD,p=0.05)
ax.set_xlabel('SALUD_std')
ax.set_ylabel('SpatialLag_SALUD_std')
plt.show()

# the map with the spots and outliers

from splot.esda import lisa_cluster
f, ax = plt.subplots(1, figsize=(12, 12))
plt.title('Spots and Outliers')
fig = lisa_cluster(lisaSALUD, 
                   MunicipiosMex,ax=ax,
                   legend_kwds={'loc': 'center left', 
                                'bbox_to_anchor': (0.7, 0.6)})
plt.show()

# significance, la significancia de cada municipio
lisaSALUD.p_sim

# quadrant: 1 HH,  2 LH,  3 LL,  4 HL
pd.Series(lisaSALUD.q).value_counts()

#Esta data la metemos al dataFrame
MunicipiosMex['ISALUD_quadrant']=[l if p <0.05 else 0 for l,p in zip(lisaSALUD.q,lisaSALUD.p_sim)  ]
MunicipiosMex['ISALUD_quadrant'].value_counts()

#Graficamos
labels = [ '0 no_sig', '1 hotSpot', '2 coldOutlier', '3 coldSpot', '4 hotOutlier']

MunicipiosMex['ISALUD_quadrant_names']=[labels[i] for i in MunicipiosMex['ISALUD_quadrant']]

MunicipiosMex['ISALUD_quadrant_names'].value_counts()

from matplotlib import colors
myColMap = colors.ListedColormap([ 'ghostwhite', 'red', 'green', 'black','orange'])



# Set up figure and ax
f, ax = plt.subplots(1, figsize=(12,12))
# Plot unique values choropleth including
# a legend and with no boundary lines

plt.title('Spots and Outliers')

MunicipiosMex.plot(column='ISALUD_quadrant_names', 
                categorical=True,
                cmap=myColMap,
                linewidth=0.1, 
                edgecolor='white',
                legend=True,
                legend_kwds={'loc': 'center left', 
                             'bbox_to_anchor': (0.7, 0.6)},
                ax=ax)
# Remove axis
ax.set_axis_off()
# Display the map
plt.show()

# The scatterplot with local info
from esda.moran import Moran_Local

# calculate Moran_Local and plot
lisaAli = Moran_Local(y=MunicipiosMex['IC_Care_en_Ali'], w=w_knn,seed=2022)
fig, ax = moran_scatterplot(lisaSALUD,p=0.05)
ax.set_xlabel('Ali_std')
ax.set_ylabel('SpatialLag_Ali_std')
plt.show()

# significance, la significancia de cada municipio
lisaAli.p_sim

# quadrant: 1 HH,  2 LH,  3 LL,  4 HL
pd.Series(lisaAli.q).value_counts()

#Esta data la metemos al dataFrame
MunicipiosMex['Ali_quadrant']=[l if p <0.05 else 0 for l,p in zip(lisaAli.q,lisaAli.p_sim)  ]
MunicipiosMex['Ali_quadrant'].value_counts()

#Graficamos
labels = [ '0 no_sig', '1 hotSpot', '2 coldOutlier', '3 coldSpot', '4 hotOutlier']

MunicipiosMex['Ali_quadrant_names']=[labels[i] for i in MunicipiosMex['Ali_quadrant']]

MunicipiosMex['Ali_quadrant_names'].value_counts()

from matplotlib import colors
myColMap = colors.ListedColormap([ 'ghostwhite', 'red', 'green', 'black','orange'])



# Set up figure and ax
f, ax = plt.subplots(1, figsize=(12,12))
# Plot unique values choropleth including
# a legend and with no boundary lines

plt.title('Spots and Outliers')

MunicipiosMex.plot(column='Ali_quadrant_names', 
                categorical=True,
                cmap=myColMap,
                linewidth=0.1, 
                edgecolor='white',
                legend=True,
                legend_kwds={'loc': 'center left', 
                             'bbox_to_anchor': (0.7, 0.6)},
                ax=ax)
# Remove axis
ax.set_axis_off()
# Display the map
plt.show()